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ABSTRACT 


A full three dimensional computational study was carried out using a finite- 
volume based solver for analyzing the performance of pin-fin based micro/meso 
scale heat exchangers with air as the working fluid. A staggered arrangement of 
cylindrical pin fins in rectangular channel geometry was used. Various 
configurations were considered consistent with a parallel experimental study 
being conducted based on a micro-wind tunnel setup. The pin/channel height 
used was 0.4 mm, and the pin diameters varied from 0.17-0.50 mm to give 
hydraulic diameters in the range of 0.13-0.78 mm. This gave volumetric area 
densities for the heat exchangers in the range of 5-15 mm2/mm3. Various heat 
exchanger configurations were simulated to determine performance 
characteristics such as the Nusselt number, friction factor, specific fluid friction 
power and Mach number in the Reynolds number regime for laminar flows. In 
addition a detailed numerical diagnosis was carried out to determine local 
behavior on the pin surfaces, end walls, etc to identify specific characteristics 
such as regions of high and low heat transfer, locations for possible shock 
formation, etc. The range of results obtained would be useful for future design of 
micro heat exchangers for use in small footprint, high heat flux dissipation 
applications like turbine blade and microelectronic systems. 
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I. INTRODUCTION 


A. BACKGROUND-MOTIVATION 

Today, great achievements and innovations are being made in the 
technological fields of microelectronics and micro-electromechanical (MEMS) 
systems. The demands for greater speed, more power, and less volume and 
mass have become more and more urgent in most of the forms and products of 
the science and technology. One of the undesirable consequences of this 
urgency is the operation in elevated temperatures. Since the systems tend to 
operate at higher energy levels, requirements are also emerging for the 
development of new devices that can remove the greater amounts of thermal 
energy and can dissipate the higher heat fluxes. The need for greater efficiencies 
and improved life cycles, which are combined with less thermal stresses, 
accelerated creep, and fatigue behaviors, is growing too. 

Two of the most concerned industries are the microprocessor- 
microelectronics and the gas turbine industries. The former has concentrated its 
efforts on dramatically reducing the size and increasing the speed of its 
attainments. This resulted in higher functional temperatures, which create a 
severe operational condition with the significant effect of limiting the devices’ life. 
So the heat removal process is critical and makes the interest in micro heat 
exchangers essential. 

The gas turbine industry faces more severe conditions in which micro¬ 
exchanger technology could be more applicable, is mainly concerned with 
increased performance which translates to a higher inlet temperature. Because 
the turbine’s reliability depends mostly on the mechanical behavior of the blades, 
their cooling due to increased performance is extremely vital. To increase the 
power, and thus the inlet temperature, the cooling techniques must be made 
more effective. 

The current methods that had been used can be divided into external and 
internal. The external method is based on air film cooling by which air is injected 
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through holes in the blade surface. The internal method is based on forced 
convection in which air is provided initially from the compressor bleed to the base 
(mounting) of the blade and then through an internal channel. A portion of the air 
can used for the film cooling, and some can be driven to a pin array, which would 
be seated at the trailing edge which is a suitable location (Ref. 1). The figure 
below presents schematically the two methods. 


Figure 1. 



Cooled turbine airfoil with pin fins [from Metzger, 1984] 


In addition to those two methods, other proposals include those that 
introduce air at the leading part of the blade, ribbed channels and impingement 
cooling techniques in order to increase the effectiveness of the heat removal 
process (Ref. 2). These two techniques dominated the research for almost two 
decades according to a literature review. The following figure shows two section 
cuts of a blade with the added characteristics, while the initial ones are kept 
parallel. 
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Figure 2. Turbine blade internal cooling passage [from Zhang et al.,1993] 


A new innovative technique, which entails implementing an array of 
staggered pins in a blade shroud and then covering it with a metal layer, was 
recently introduced. The air after its passage through the micro-exchanger 
passages can be driven outwards for external film cooling. Two manufacturing 
attempts have been reviewed; one by electro-chemical deposition (Ref. 3), the 
other one by LIGA process (Ref. 4). 



Figure 3. SEM of the heat exchanger [from Marques, 2002] 
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Shroud 




Figure 4. Cross section of heat exchanger in stainless steel tube [from 

Marques, 2002] 


The above figures present the configuration of the exchanger mounted on 
a model blade and the operation during which air flows through the array and 
cools the hot blade. 

A definition sketch of the chosen micro-heat exchanger array geometry 
consisting of cylindrical pins is shown in figure 5. In the geometry, X stands for 
the streamwise pitch between two pins, S is the span-wise pitch, H is the pin 
length/height and D is the pin diameter. This nomenclature is consistent with the 
literature. 

Flow 

- > 


Flow 

-> 

Figure 5. Fleat exchanger scheme with the common symbolism [from Choo, 

2003] 
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B. LITERATURE REVIEW 

Most of the work in this area so far has been empirical and mainly in the 
macro-scale. Several researchers tried to understand the performance of 
different array configurations, in several conditions, by an examination of factors 
like the Nusselt number, the heat transfer coefficient and the friction factor. The 
latter is as significant as the other heat transfer characteristics, since in order to 
evaluate a configuration, the pressure drop has to be considered. The usual goal 
in these studies is to optimize the array to provide the greatest heat exchange 
rate with the least expended work for any flow condition. 

One of the early contributions came from Van Fossen (Ref. 5), who 
studied the effect of the presence of pins with different H/D values, for a range of 
Reynolds number from 3,000 - 60,000. He concluded that the presence of the 
array returned higher heat transfer coefficient values compared to an empty duct 
and they got higher as the H/D was increased. His results correlated better flows 
over a Reynolds number of 6,000. 

Metzger et al. (Ref. 1, 6) developed accurate correlations for 
configurations in which the stream-wise distance ratio was varied while the 
others were kept constant. Those correlations for the Nusselt number were 
validated by Hamilton (Ref. 7) with a numerical model. Also Metzger was one of 
the few who tried to examine the effects of using other types of pins, like oblong 
shaped, which was also the research subject of Arora (Ref. 8). Both Metzger and 
Arora concluded that the friction factor was significantly lower, but the heat 
transfer rates did not approach the cylindrical pin values. 

Chen et al. (Ref. 9) later examined more aerodynamically shaped pins like 
tear drop shaped fins. They concluded that these pins offer higher performances, 
and better heat transfer rates, with almost half the friction factor of the round 
pins. Those findings were corroborated by Hamilton, who found that the optimum 
airfoil shaped pin array outperforms a similar cylindrical pin array, which needs 
triple the specific energy loss to produce analogous heat transfer rates. 
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There has been considerable deliberation among the researchers with 
regard to the heat fluxes from the pins compared to the endwalls. Chyu et al. 
(Ref. 10) show comparable fluxes between these two parts of the exchanger; 
unlike Metzger and Van Fossen, who showed that the pins generated 50 percent 
higher heat transfer rates than the end-walls. On the other hand, Al Dabagh et al. 
(Ref. 11), found that end-walls show heat transfer rates 35 percent higher than 
the staggered pins. Finally, Chyu et al. (Ref. 12) in alter study found that the pin 
surfaces offer around 20 percent more in heat transfer than the endwalls. 
Furthermore their results for the Nusselt number are very close to those reported 
by Metzger (Ref. 6) 

Also very interesting is the approach of Tahat et al. (Ref. 13, 14), who tried 
to correlate the Nusselt number with the Reynolds number, including the 
geometry configuration of the tested arrays by using the ratios of the stream and 
the span-wise distance of the pins over the total length and width of the array. 
They used a specific test bed with a varying ceiling height to test different pin 
bank arrangements to find the optimum spacing for the pins. 

In 1997, Li et al. (Ref. 15) conducted an experimental work with elliptical 
pin arrays and compared their data with Metzger’s findings. They used the same 
dimension pins in two different arrays, but with equal stream and span-wise 
distance between them, for a range of the Reynolds number 800 - 9,000. They 
found better heat transfer rates with considerable lower pressure losses. The 
more packed his configuration was, the greater the performance with resulting 
higher Nusselt numbers. 

One of the most recent, and probably the first in micro dimensions 
experimental work is by Marques et al (Ref. 16), who tested the configuration of 
X/D=S/D=2.5 and Fl/D=1 for Reynolds number from 4,500 to 19,200. The same 
configuration was also tested by Metzger and Chyu in the past, but in the macro¬ 
scale. Marques et al. developed a correlation for the Nusselt number for the 
range of Reynolds noted above, which gives more optimistic results than the 
other two researchers. A possible explanation is that their theoretical model uses 
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a smaller heat exchange area in the array than the exchanger footprint area. 
They then developed a theoretical model that can compute the performance of 
such a heat exchanger when it is used under a shroud in a turbine blade to 
cooling it. 

Finally, at the Naval Postgraduate School, Monterey, California, in 2002- 
04 an experimental apparatus was built for testing a set of compact heat 
exchangers by Roussakies (Ref. 17) and Bosserman (Ref. 18). The arrays are 
the same as those that will be studied in this work. Because some design issues 
of the testing bed facility arose, the results are not as reliable as had been 
expected which a good reason for also checking them against a numerical 
model. 

On the other hand, little numerical modeling work had been done until 
recently in comparison with empirical. Recent advances in CFD tools have made 
this a very promising avenue of research. Two important factors in the push 
toward that direction are the great manufacturing and testing difficulties and the 
high cost that these studies demand. 

Donahoo et al. (Ref. 19) developed a two dimensional model of staggered 
pin arrays, by using a general purpose viscous solver, in order to optimize their 
configurations in axial pitch (X/D) spacing for turbulent Reynolds numbers. They 
found that the maximum heat transfer occurred between the fourth and sixth row 
for the tested arrays. Finally they recognized that a 3 D model is a future 
necessity since it was not possible to capture the pin-endwall interactions with 
the 2 D model and compare their participation to the whole heat transfer of the 
array or to study the effect of the FI/D ratio. 

Flamilton (Ref. 7) used the finite element commercial code “ANSYS” to 
develop a numerical 3 D model to study several configurations of staggered pin- 
fin heat exchangers, but still in the macro-scale. Initially, he validated his model 
with prior experimental results from Arora et al. and Metzger et al. Fie then 
proceeded to examine of the performance of arrays with cylindrical pins. Fie 
developed correlations for the Nusselt number for all his configurations, and also 
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found that pin surfaces had from 50-200 percent higher heat transfer coefficients 
than the walls. Achieving small increases of it requires disproportional increases 
in friction power losses. He concluded that variations in axial pitch did not have 
as great an effect in performance as those the in H/D ratio. He also studied 
airfoil-shaped pin arrays, from which he concluded that they are better 
performers than corresponding dimensional cylindrical pin arrays since, to attain 
comparable heat flux removal, they require only a third of the energy as 
compared to the latter. 

Choo (Ref. 20) carried out a laminar flow numerical analysis in 2003. He 
simulated exchangers performing in micro-scale in the laminar regime with 
Reynolds numbers in the range of 100-1,000 to study the effect of the 
geometrical ratios of the exchangers (X/D, S/D, H/D). He stated that for his 
geometries and simulations the maximum heat transfer occurred around the 
second row of the models. Also, after taking into account all three effects, he 
suggested that the optimum configuration of his matrix was X/D=1.25 and 
S/D=H/D=3. 

To minimize the amount of time and cost that are needed in an empirical 
project, use of a computational study for pre-evaluation and optimization of 
several designs is absolutely essential. Another characteristic worth mentioning 
is the tremendous flexibility that a computational study can offer as compared to 
empirical methods. Finally, the attributes of the theory that were developed for 
the macro-scale studies can also be directly applied to micro-scale, since they 
are based on the same characteristic qualitative factors. 

C. OBJECTIVES 

The main objective of this study is to determine the performance and 
characteristics of flow and heat transfer behavior of a set of micro heat 
exchangers of specific configurations. A 3-D model will be developed that will be 
validated by theoretical and historical experimental results in order to use it to 
examine the behavior of the different staggered cylindrical pin arrays in a laminar 
compressible flow regime. 


8 



The parameters that will be inspected are the Nusselt number, the heat 
transfer coefficient, the Mach number, the array pressure drop, the friction factor, 
and the specific fluid friction power for different Reynolds numbers The 
theoretical definitions for the above parameters will be provided in the following 
chapter. The results will be used to define suitable regions-conditions of 
operation for every heat exchanger’s geometry and will be compared to findings 
of experimental tests already completed. 


D. METHOLOGY-ORGANIZATION 

Several three-dimensional models were constructed according to our test 
matrix, using the commercial, computational fluid dynamics package CFD-ACE+. 
The construction was carried out with the preprocessor of the package CFD- 
GEOM. Then the models were transferred to the processor CFD-ACEU for 
executing the simulation. Finally, the postprocessor of the package CFD-VIEW 
was used to study, analyze, and manipulate the results. The details of the 
modeling will be discussed in the next chapter. In the third chapter we provided 
the validation and corroboration with the theoretical and historical results and a 
more extensive heat exchange examination of the observations. In the fourth 
chapter the results of the simulations of the tested configurations are discussed. 
The final chapter draws certain conclusions and provides some 
recommendations. 
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II. FINITE VOLUME MODELING 


A. MODELING GENERAL DESCRIPTION 

The CFD-ACE+ package solver was used in this study. It is a solver that 
uses a finite volume pressure correction based method to solve the conservation 
equations of mass, momentum and energy. The software provides the user with 
three graphical user interfaces (GUIs): GEOM for the geometries-models 
creation, the ACEU that is the actual solver of the discretized equations, and the 
postprocessor VIEW for viewing and post-processing the results. CFD-ACE+ has 
been very popular for this kind of problem, such as fluid flow fields in conjugation 
with heat transfer, because it provides reliable results for a variety of problems. 
This independence makes this solver more powerful since it provides very similar 
results for several grid refinements within good accuracy. 


B. TEST MATRIX 

The test matrix of this study consists of ten different configurations. The 
first eight configurations and the tenth have been tested here at the Naval 
Postgraduate School by testing Bosserman and Roussakies; but because of 
some construction constraints and limitations of the utilized bed, their data is still 
preliminary and is being verified. All of the configurations have the same height 
due to the specific height of the flow channel inside the designed apparatus 
(Fl=405pm). The first four have a pin diameter of 500pm (D=500pm), while the 
next set of four has a pin diameter of 166.7pm (D=166.7pm). The tenth one has 
a pin diameter of 151.5pm; it was Summers (Ref. 21) optimum laminar 
configuration in the macro-scale level. The two sets of four vary under the same 
ratios for the stream-wise and the span-wise pin spacing. This can be observed 
in the following table, which presents the tested matrix, with some geometric 
characteristics that will be explained thoroughly later in this chapter. 
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Table 1. Test matrix for the heat exchanger configurations 


# 

X/D 

S/D 

H/D 

D (pm) 

oH 

Dh(pm) 

Vopen(ni ) 

Awetted(ni ) 

HX#1 

3 

3 

0.81 

500 

2.11 

701 

8.32E-10 

4.74E-06 

HX#2 

1.25 

3 

0.81 

500 

2.26 

567 

3.00E-10 

2.12E-06 

HX#3 

1.25 

1.25 

0.81 

500 

2.62 

307 

7.87E-11 

1.02E-06 

HX#4 

3 

1.25 

0.81 

500 

2.26 

567 

3.00E-10 

2.12E-06 

HX#5 

3 

3 

2.43 

167 

2.67 

553 

9.24E-11 

6.68E-07 

HX#6 

1.25 

3 

2.43 

167 

3.62 

354 

3.34E-11 

3.77E-07 

HX#7 

1.25 

1.25 

2.43 

167 

5.88 

137 

8.74E-12 

2.55E-07 

HX#8 

3 

1.25 

2.43 

167 

3.62 

354 

3.34E-11 

3.77E-07 

HX#9 

2.5 

2.5 

1 

500 

2.25 

777 

6.83E-10 

3.52E-06 

HX#10 

5 

5 

2.67 

152 

2.27 

690 

2.25E-10 

1.30E-06 

The 

ninth 

configuration 

is the 

one 

used most by the 

researchers. 

especially 

in the 

macro-scale 

and it is 

1 the 

one that 

is going to 

be used for 


validations and corroborations, along with the seventh, which has the largest 
magnification factor (oH) and is the “densest” of all of them. The area density of 
an array unit cell is defined as the ratio of the wetted (heat transfer) surface area 
over the total volume. The magnification factor that is provided in the table is 
defined as the product of the area density and the pin height. For an exchanger 
without pins that is a rectangular cross sectional duct, the oH becomes equal to 
2. All the above table values have been computed based on a unit cell. 


Table 2. Unit cell dimensions for the tested heat exchangers 


# 

X/D 

S/D 

H/D 

D (pm) 

S (pm) 

X (pm) 

H(pm) 

HX#1 

3 

3 

0.81 

500 

1500 

1500 

405 

HX#2 

1.25 

3 

0.81 

500 

1500 

625 

405 

HX#3 

1.25 

1.25 

0.81 

500 

625 

625 

405 

HX#4 

3 

1.25 

0.81 

500 

625 

1500 

405 

HX#5 

3 

3 

2.43 

167 

500 

500 

405 


12 



HX#6 

1.25 

3 

2.43 

167 

500 

208 

405 

HX#7 

1.25 

1.25 

2.43 

167 

208 

208 

405 

HX#8 

3 

1.25 

2.43 

167 

208 

500 

405 

HX#9 

2.5 

2.5 

1 

500 

1250 

1250 

500 

HX#10 

5 

5 

2.67 

152 

758 

758 

405 


The dimensions of the unit cells for all the configurations are shown in the 
Table 2. The heat exchangers have been divided into two groups and they will be 
examined according to their pin diameter. 


C. BOUNDARY CONDITIONS 

Several assumptions had to be made for the purpose of reducing the 
computational requirements. Because the geometry and scale of the micro heat 
exchanger it usually contains a large number of pins, and the contribution of pins 
as compared to the endwalls in heat transfer will be sufficiently greater. Micro 
heat exchangers have a very large number of pins in the span-wise direction, 
and the contribution of the side walls in the heat exchange process could be 
relatively neglected. Instead Hamilton found that it is sufficient to consider only 2 
spanwise rows and introduce symmetry sidewall boundary conditions to account 
for the ignored rows. A symmetric boundary was also assumed at the Hx 
midplane in order to simulate only half of the geometry. The symmetric planes 
have the zero cross-over velocity components, and from the thermal perspective 
is are adiabatic in nature. Furthermore, the wall boundary condition was used for 
the pin surfaces and the lower endwall where the no slip condition was applied. 

From the heat transfer point of view, the pins and the lower endwall were 
assumed to be isothermal. A justification calculation for the isothermal character 
of the pins can be found in Hamilton (Ref. 7). An illustration of the model is 
presented in figure 6. 
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Heated pin 
Surface 


Isothermal 

Heating 


Figure 6. Schematic of the numerical model 


To validate the models with the theoretical and historical experimental 
results, an entry and an exit straight duct section was used. The entry duct is 
used to provide a fully developed velocity profile at the flow entrance to the 
exchanger; while the exit duct will provide well mixed conditions for the flow so as 
to make the exit temperature calculations more effective and robust. 
Subsequently, we used the concept of the entry and exit ducts in all the array 
configurations in order to be consistent. In both ducts we applied the same 
boundary conditions for the surrounding planes. The sidewalls and the upper 
plane were assumed to be adiabatic symmetry planes, while the bottom was 
assumed to be an adiabatic wall. Furthermore, we assumed an inlet with a fixed 
velocity that fixes the mass flow rate and thus the Reynolds number of the 
simulation and a fixed pressure outlet boundary. A representation of the above 
conditions is provided in the figure 7. 
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Heated pin surface 
at T=320K 


Heated endwall 
at T=320K 


Air inlet 
T=300K 



Adiabatic exit endwall 


Entry duct 


Adiabatic Symmetry Condition, Uy=Uz=0 


Adiabatic entry endwall 


Figure 7. Boundary conditions representation 


D. GRID GENERATION 

Since we decided to use the “GUI” of the GEOM and not the PYTHON 
scripting feature, we had to create several geometrical models independently. 
First, the 2-D model of a half unit cell had to be constructed, since we intended to 
use an unstructured triangular meshing. Thus the creation of a closed loop was 
unavoidable. The next step was to duplicate that half unit cell at the proper 
position, according to every array configuration. We did a very careful meshing 
that was finer near the pin’s perimeter and coarser towards the internal part, to 
capture the large heat and flow gradients near the walls. In most cases we used 
a 30 degrees curvature criterion and a growing factor from 1.05 to 1.10. Then we 
linearly extruded the whole geometry toward the third direction by using a power 
law with the purpose again to create layers of triangular prisms that will be finer 
near the heated endwall where the gradients are also expected to be higher. We 
carefully examined all the 3-D centroid / face angles to ensure that were at least 
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25 degrees and also that no negative volume cells had been created. Figure 8 
provides a fingerprint illustration of the grid around a unit cell of an array. 



Figure 8. Sample unit cell’s meshing in 2-D 


Another issue was the entry length in the ducts. Since we had to construct 
several configurations and run them in different Reynolds numbers we had to 
figure out a way to reduce the model generation work and eliminate the need for 
different ducts for every simulation. This was achieved by using a feature of the 
solver ACE where we can simply input the velocity profile for the inlet boundary 
to ensure fully developed conditions. Therefore, we used a constant length duct 
for all the models that were generated, with lengths of three and five hydraulic 
diameters for the entry and exit sections respectively. Since the height of all the 
models is the same, the duct hydraulic diameter for all the models is constant 
and equal to twice the height for an infinite span-wise model. 
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Because we used structured grids for meshing the two ducts, the whole 
model had a hybrid character, and we tried not to exceed a block aspect ratio of 
five in the cell extrusion process. This would serve the connection regions 
between the two grids for having similar cell dimensions to ensure a smooth 
transition from one volume to the next. 



Figure 9. Sample model’s grid 


Figure 8 illustrates a sample mesh of a model. After the generation of any 
geometry we had to scale it in micrometers in order to import that into the solver 
GUI. 

E. SOLUTION METHOD - SOLVER CONTROLS 

As mentioned previously the simulations were executed in the laminar 
regime. Hamilton (Ref. 7) in his dissertation showed that the transition typically 
occurs at a Reynolds number 1000-1500. Within that region, he achieved an 
overlap in a figure of the laminar and turbulent models that he used. Zhukauskas 
(Ref. 23) makes a similar statement in his research review. The Reynolds 
number on which our runs are based extended from 100 to 1,000. 
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In the CFD-ACE+ solver the flow and heat transfer modules were 
selected in order to obtain the flow (velocity and pressure) and temperature 
distributions for the fluid volume elements. To formulate an equation for the 
unknown variable pressure, the continuity equation is solved by adopting the 
SIMPLEC algorithm in CFD - ACE to start the iteration process. For more 
information, the reader is referred to the user manual (Ref. 24) or the book by 
Versteeg and Malalasekera (Ref. 25). 

In regard to the volume conditions, because we are going to solve the 
compressible form of the numerical models, the density is defined from the ideal 
gas law, which is the sixth equation in the governing equation set. Furthermore, 
the air viscosity was obtained by using Sutherland’s law for. Also, the air 
pressure coefficient Cp was held constant, since its variation in the region of 300 
- 320 K is negligible. For the thermal conductivity, which varies more than the 
pressure coefficient, we used the constant Prandtl number law which is a 
variation of k according to the definition of the Prandtl number which is held 
constant at a value of 0.707. 


Pr = 



^k = 



( 1 ) 


The next step in the solver tabs is to input the boundary conditions. The 
velocities for the wall boundaries, pins and endwall were kept at zero and the 
temperature was set to 320 K. The inlet temperature was set to 300 K, and for 
the flow part we input a fully developed X-velocity profile in the Z direction. The 
profile for every Reynolds number and for every array was calculated by using an 
ECXEL spreadsheet; which was then utilized by the solver to fix the mass flow 
rate which subsequently derives the Reynolds number for this simulation. 

For the flow initial conditions we input a velocity that was based on the 
duct cross sectional area and was derived by dividing the specific mass flow rate 
with by duct area and film density, while the initial temperature was set to 300 K. 

Another crucial control is the scheme that we had to choose in order to 


discretize the advection term, especially when this term dominates other terms in 
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the governing equations. The monotone streamline upwind (MSU) scheme is 
generally quite robust and produces stable diagonally dominant matrixes, but is 
first order accurate. It is also the scheme recommended in the ACE+ manual, 
especially for complex geometries. 

Trial runs were carried out using other second order schemes like the 
central, or the upwind limiter, always blended with the MSU to prevent 
divergence. The same scheme is also suggested by the solver manual 
documentation. We tried to overcome the choice of the first order scheme by 
increasing the convergence criteria for the residuals average and by using more 
refined models for the simulation. As an example of the deviation between the 
schemes for HX#7, and for the case of Re=100, the percentage difference 
between the upwind and the central, as far as it concerns the heat addition to the 
fluid was 0.2 percent and between the upwind limiter was 0.3 percent. 

As far as it concerns the relaxation criteria, in most cases the solver 
default values worked well, except for the higher Reynolds number simulations, 
near 1,000, where we had to modify them a little by increasing the inertial 
relaxation values up to 0.3 - 0.4 and by decreasing the linear ones down to 0.9 in 
order to provide more stability to the solution. This slowed down the convergence 
process but made it achievable. 

F. NUMBER OF PIN ROWS 

Another issue that we had to resolve was the number of pin rows that we 
were going to use for every array that we were going to model. Because of the 
great variation of the configurations, a great variation in their performances was 
expected. For every configuration, and especially for the lower Reynolds number 
simulations, the bulk outlet temperature was expected to be higher. It began to 
reduce as the Reynolds number was increased, which agrees with the theory, in 
Kays and Crawford “Convective Heat and Mass Transfer” (Ref. 26). So for every 
array we had to use only so many rows of pins that would yield an outlet bulk 
temperature no more than 319 K to 319.5 K, or 0.5 to 1 K lower than the 
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prescribed wall and pins temperature of 320 K, in order to avoid errors due to 
loss of resolution in the numerical calculations. 

Therefore, for every configuration we generated a geometry without ducts 
and with ten rows of pins, and we ran that for Reynolds of 100. Then by utilizing 
CFD - VIEW, we constructed a carpet plot of the temperature for the upper 
symmetry plane, illustrated in figure 10, to get an initial rough estimate in which 
row we achieved a temperature around 319 K. 



Figure 10. Temperature profile for two pin arrays with ten rows 


Then for higher accuracy we went to this row of fluid blocks, made a 
stream wise direction (X) cut and integrated the results on that plane. This was 
done, to find an accurate bulk temperature at that exact X coordinate. By moving 
the cut in the selected row blocks we could check the temperature continuously 
in every position. 

Then we went back to GEOM to generate a final version of the model with 
the desired number of pin rows and with the entry and exit duct as it were 
described above for the specific runs. 
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Table 3. Number of pin rows used per configuration 


# 

X/D 

S/D 

H/D 

Nx 

Ns 

L(m) 

W(m) 

HX#1 

3 

3 

0.81 

6 

2 

0.0090 

0.0030 

HX#2 

1.25 

3 

0.81 

10 

2 

0.0063 

0.0030 

HX#3 

1.25 

1.25 

0.81 

4 

2 

0.0025 

0.0013 

HX#4 

3 

1.25 

0.81 

4 

2 

0.0060 

0.0013 

HX#5 

3 

3 

2.43 

8 

2 

0.0040 

0.0010 

HX#6 

1.25 

3 

2.43 

10 

2 

0.0021 

0.0010 

HX#7 

1.25 

1.25 

2.43 

4 

2 

0.0008 

0.0004 

HX#8 

3 

1.25 

2.43 

4 

2 

0.0020 

0.0004 

HX#9 

2.5 

2.5 

1 

6 

2 

0.0075 

0.0025 

HX#10 

5 

5 

2.67 

8 

2 

0.0061 

0.0015 


Table 3 provides the specific number of rows that every heat exchanger 
configuration used. All the simulations used two unit cells in the span-wise 
direction. In the above table the dimensions of the configurations for all the heat 
exchangers are also presented. 


G. THEORETICAL APPROACH 

1. Characteristic Length 

The most appropriate characteristic length that is the hydraulic diameter of 
this Hx array is obtained from flow theory in porous media, and is defined as 
follows in terms of the open volume that is available for the fluid flow in the array 
and the wetted heat transfer surface area: 

4*V 

( 2 ) 

^wetted 

It is the most truly dimensionless quantity that represents consistently the 
tortuous flow paths inside the exchanger. It is the same characteristic length that 
was proposed by Van Fossen (Ref. 5) and that is used by many more recent 
experimental and numerical researchers. The hydraulic diameter remains the 
same whether its calculation are made based on the unit cell or based on the 
whole array. 
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2. Entry - Exit Length Considerations 

According to White (Ref. 22), for every Reynolds number the relation that 
gives the entrance length for a fully developed flow is the following: 

L3«0.06*D,*Re,, (3) 

Because it was not possible to construct a different entry duct for every 
simulation and because the existence of ducts usually smoothes the flow 
conditions, we decided to impose a fully developed velocity profile and to use the 
constant length specified above for all the flow conditions. 

For the exit duct length, which is more crucial because it is necessary for 
ensuring a well mixed flow that is starting to redevelop again and because it will 
help in a more accurate integration for computing the predicted outlet bulk 
temperature, we used a little longer section than the entry one, equal to five 
hydraulic diameters. 

3. Inlet Velocity Considerations 

According to White (Ref. 22), the relation that gives the velocity profile for 
a Poiseuille flow between parallel plates is 

( 4 ) 

where Umax is the maximum velocity at the centerline of the duct, which is 
equal to 

U™„=|V, (5) 

where Vin is the average inlet duct velocity. Consequently by using the 
equation (4) we created a six - point profile along the z-axis, which then was 
used at the definition of the inlet boundary condition in the solver. 

4. Reynolds Number 

The Reynolds number that we are going to use is based at the previously 
defined hydraulic diameter and in a characteristic velocity, the average array 
velocity 
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( 6 ) 


Re = 


pLJD, 

M 


where 



(7) 


and A is the average array cross sectional area, which is equal to 

- V 

A = (8) 

Also, according to the above relations, the Reynolds number can be 
expressed associated to the mass flow rate. This relation will be used after every 
simulation to calculate the exact Reynolds number from the “mass flow 
summary” from the output file of ACE. 


Re = 


pA 


(9) 


H. COMPARED PARAMETERS ANALYSIS 

1. Nusselt Number 

The Nusselt number is a dimensionless number that provides the ratio of 
the convective heat transfer over a surface that would occur by fluid motion to the 
corresponding conductive heat transfer and is a measure of the heat transfer 
ability. In our case more meaning makes the average Nusselt number which is 
given by the following relation: 

= ( 10 ) 

where hg^^gy is the average array heat transfer coefficient and can be 

calculated by using the thermal energy difference between the flow inlet and the 
outlet, which is provided by the heat transfer summary in the output file of ACE at 
the end of the simulation, and by using the log mean temperature difference as 
below: 
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AQ 


( 11 ) 


h 


array 


A AT 

'^wetted'-' 'im 


For computing the log mean temperature difference we will use the 
definition of the outlet bulk fluid temperature, which, according to Kays and 
Crawford in “Convective Heat and Mass Transfer” (Ref. 26) is provided from 


rfiT 


bulk,out 



( 12 ) 


In the above definition the Cp is not included because we had assumed 
already that it is constant and does not influence the relation. The outlet bulk fluid 
temperature was calculated in the VIEW program by integrating at the outlet of 
the adiabatic exit duct with respect to the duct cross sectional area by using the 
expression calculator feature. It is the temperature that characterizes the average 
thermal energy of the flow, since it is based at the mass flow rate. 

The log mean temperature difference is given by the formula below; 


A-p _ ("l”wall "l"bulk,in )-{\ 
fT _T 


where the inlet bulk fluid temperature was constant and equal to 300K. 

2. Effective Heat Transfer Coefficient 

The effective heat transfer coefficient is the corresponding heat transfer 
coefficient in an empty (without pins) open channel had to provide in order to 
transfer the same amount of thermal energy like the one with the pins. 

u u aH .... 

heff =harray^ (^4) 

It is an average quantity of the whole array since it is based at the average 
array heat transfer coefficient, and it becomes equal to it for the open channel 
because aH=2. 
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3. Friction Factor 

The dimensionless friction factor of the pin array that was used was 
defined as follows: 

f ^ AParrayDh 

-pU\ 

2 

where AParray was the total pressure difference between the inlet and 
outlet sections of the array. It was computed from the VIEW post processor by 
integrating at the sections with the expression calculator tool in order to find an 
average value of the inlet and outlet section. L is the total stream-wise length of 
the array. Also, the density that was used was the one based on the standard 
atmospheric conditions but calculated in the film temperature, in order to have 
the same measure of comparison between the different heat exchanger 
configurations. 

4. Specific Fiuid Friction Power 

Kays and London (Ref. 27) define a quantity that can express the 
performance of a heat exchanger, the specific fluid friction power, which 
characterizes the work that must be expended to overcome the viscous effects 
per unit wetted area. 

E= (16) 

Pfilm^ wetted 

5. Effectiveness - NTU Reiations 

According to Incropera and De Witt in “Introduction to Heat Transfer” (Ref. 
28) the effectiveness of a heat exchanger is defined as the ratio of the actual 
heat transfer rate to the maximum possible heat transfer rate and, in our case, is 
defined as: 


£ = 


"^bulk.out "^bulk.in 
"^wall "^bulk.in 


(17) 
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The same source shows that, for all the exchangers, the effectiveness 
associated with the number of transfer units (NTU) is given by the following 
relation: 


£=1-exp(-NTU) (18) 
The NTU generally is defined as 

l lA 

NTU = — (19) 


but in our case as 


1 

UA 


array 


A 


wetted 


(20) and C,,=mCp (21) 


6. Maximum Mach Number 

The maximum Mach number was extracted with the post processor by 
examining the fluid volume inside the array and by again using the expression 
calculator and the max built in function. In order to define suitable operating 
regions for all the exchangers, the maximum Mach number will be plotted versus 
the corresponding Reynolds number for comparison. 
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III. VALIDATION - CORABORATION AND TRANSITION 

EXAMINATION 


A. GRID INDEPENDENCE STUDY 

It was anticipated that the simulations of these micro-scale heat 
exchangers would yield greater performance values when compared to their 
macro-scale counterparts. This meant that much higher gradients of the 
concerning parameters were expected to exist within the solutions. Therefore, 
the need for a suitable grid arose, and this need determined how fine the grid 
should be. It was a significant effort to obtain valid (grid independent) results. The 
following section indicates the procedures taken to achieve grid-independent 
results. 

Two prismatic elements were considered for the internal part of the 
exchangers: a quadrilateral element and a triangular element. Following the 
element selection, a specific measure was identified to properly quantify the 
quality of such grids. The most suitable measure for this task was the pin 
perimeter discretization. Furthermore, to make a consistent grid comparison, the 
grid growing factors toward the planar and orthogonal directions to the pins’ 
cross sectional area were kept constant. Analogous refinement was also 
conducted within the entry and exit sections of the heat exchangers. 

The grid independent study was performed with the X/D=1.5, S/D=1.5 and 
Fl/D=1 heat exchanger configuration with a pin diameter of 100pm. This 
configuration was selected due to its high area amplification factor, which 
indicates a very high pin density. It was anticipated that the solutions would have 
high gradients in all concerning parameters, which gives the solution a very 
strong dependence on grid refinement. Thus, to achieve grid independence the 
computed heat transfer flux had to reach less than two percent variation between 
subsequent grid refinements. 


27 



Table 4. 


Grid independence for model with triangular elements 


Parameter 

Grid 1 

Grid 2 

Grid 3 

Perimeter division 

50 

70 

90 

# of cells 

88143 

141076 

201165 

# of nodes 

137380 

159516 

234304 

AQ (W) 

0.00271 

0.0027 

0.00268 

Error (%) for AQ 

1.12 

0.37 

0.00 

T bulk out (K) 

305.49 

305.46 

305.44 


Table 5. Grid independence for model with quadrilateral elements 


Parameter 
Perimeter division 

# of cells 

# of nodes 
AQ (W) 

Error (%) for AQ 
T bulk out (K) 


Grid 1 

Grid 2 

50 

70 

87450 

139802 

72030 

118800 

0.00257 

0.0026 

4.10 

2.24 

305.21 

305.32 


Grid 3 

Grid 4 

90 

110 

198660 

278018 

172284 

245504 

0.00265 

0.00268 

1.12 

0.00 

305.37 

305.43 


Tables 4 and 5 show the pin perimeter discretization characteristics and 
the significant data extracted from the simulations. The triangular elements were 
selected because they proved to be more robust and because the solutions with 
the triangular elements required less grid refinement to obtain stable outlet bulk 
temperatures. 


28 



AQ for testing T=306K for different grids 



Figure 11. Comparison of heat transfer vs grid quality for the examined 

models 

Figure 11 compares the results from simulations performed with both, the 
triangular and the quadrilateral element models versus their grid quality (number 
of discretization points within one pin perimenter.) These tests were performed 
for the case of Reynolds number of 100. As shown, both models converged to 
the same solution of 0.00268 W. Flowever, the triangular modeling show results 
that overestimates the heat transfer when using a grid too coarse to properly 
capture the high gradients. In contrast, the quadrilateral elements underestimate 
the heat transfer. Therefore, selecting the triangular elements was simple, 
specially noting that grid independent results were obtained with 20% less 
discretization points around the pin perimenter. Furthermore, the deviation from 
the grid independent result was 300% larger with the quadrilateral elements. 

Figure 12, plots the models’ performance versus the time required for the 
simulation execution. Selecting the number of divisions per pin perimeter was 
also important. From both plots (figures 11 and 12,) the choice of 70 points 
provided the best trade-off between computational time and accuracy of the 
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results. At 70 points per pin diameter, the results predicted a heat transfer rate of 
only 0.4% different than the finest configuration (using 110 points). 



Figure 12. Comparison of the outlet temperature achieved vs. computational 

time 


B. THEORETICAL VALIDATIONS 

Numerous inspections were performed in order to validate the results 

generated by CFD - ACE. The contour plots of temperature, velocity and 

pressure were examined carefully to verify that they satisfy the boundary 

conditions. The velocity profiles were inspected for a fully developed condition 

before the entrance into the array. In addition, the output file generated upon 

completion of every simulation was carefully examined and analyzed. For 

example, the mass flow rate computations were examined for an exact 

agreement with the inlet and outlet boundaries. This was verified by assuring that 

the imbalance was smaller than 1% of the total value. The heat transfer summary 

was also carefully inspected. The energy flow into and out of the model at the 

corresponding boundaries was examined to be equal to the total amount added 
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from the pins and the endwall. Also, the heat imbalance was examined to check 
if it was an order of magnitude smaller than the net energy addition. However, 
because the magnitude of the non-dimensional parameters is extremely small 
(due to the heat exchangers’ dimensions,) computation errors were introduced. 

1. Energy Balance 

This was examined for agreement between the bulk outlet temperature 
during the post-processing of the results and the temperature that was computed 
from a simple energy balance according to the following relation: 


TDut Tn 


+ ^ ( 22 ) 
mC„ 


The AQ was computed from the heat transfer summary from the output file 
runoff each simulation. If the two temperatures had more than a 0.1% difference, 
the simulation was repeated with a finer grid or with an increased number of 
iterations. In most cases, we did not have any difficulty, except for a few 
simulations at a Reynolds number of 1,000 and for some of the configurations 
with the highest area magnification factors. 

2. Effectiveness 

The independence of heat exchanger effectiveness was tested by varying 
the prescribed boundary temperature conditions. Thus, during the grid 
independence tests, three different boundary temperatures, 306 K, 310 K, and 
315 K were used with all the grid refinements. The results, shown in table 6, were 
extremely close with values within 0.1% from each other. 


Table 6. Effectiveness independence of boundaries’ temperatures 


Model 

1st (50 points) 
2nd (70 points) 
3rd (90 points) 
4‘^ (110 points) 


Twall,1~306 

0.868 

0.887 

0.896 

0.905 


Twall,2~31 

0 

0.868 

0.887 

0.897 

0.905 


Twall,3~ 315 

0.869 

0.888 

0.897 

0.905 


31 



For grater assurance, heat exchanger No. 7 was used using two very 
distinct boundary temperatures, 306 K and 400 K. The effectiveness results were 
again within 0.1% of each other. 

3. Effectiveness Based on NTU 

Effectiveness is associated to the number of transfer units (as shown in 
equation 18). Since this relation is valid for every heat exchanger, it can be used 
to further validate the computational results. In order to thoroughly test the 
computational model against the trend of this relation, heat exchanger No. 7 was 
used. Two different models were created, one with two rows of pins and the 
second with four rows of pins.. The computed NTU results are plotted in the 
purple squares of figure 13. 



Figure 13. Effectiveness based on NTU for FIX#7 
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The black line in Fig. 13 plots the computed NTU values using Eq. 18, the 
results show how close the CFD models approach the theory. 

4. Compressibility Effects 

Our initial ideas of simulating the problem by using a constant density for 
the models, similar to previous numerical works, collapsed very early since some 
of the pin-fin arrays are highly packed. It was anticipated that these arrays would 
require high differential pressures across the array, which breaks down the 
assumption of constant density. 

Fleat exchanger configuration No. 7 was used to simulate various mass 
flow rates with and without variable density configurations. The tests were 
performed at Re = 100 and 1,000. The variable density model used the ideal gas 
law to compute the appropriate density values. Figure 14 shows these results, 
which show significant difference in their friction factors, up to 50% difference at 
Re = 1,000. 


Effect to "f due to density variation vs Reynolds number 



Re 


Figure 14. Comparison between constant and variable density for HX#7 


Consequentially, variable density was incorporated in all subsequent tests. 
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C. LAMINAR TO TURBULENT TRANSITION 

An examination of the transition region for the micro heat exchangers was 
also necessary. For this case, heat exchanger configuration No.9 was used, 
because the results can be corroborated with existing experimental data. 

1. Turbulence Modeling 

Within our computational tool, CFD-ACE, various turbulence models can 
be used. The more often used, the “kt-£t” model, does not offer very accurate 
results during transition because it only considers the overall effect of turbulence 
in the mean flow field. Flowever, when the interest occurs with parameters such 
as the heat transfer coefficient and the friction factor during laminar-to-turbulent 
transition, it is suggested to use a Low Reynolds number model. These models 
permit the integration of momentum and “kt-£t” equations all the way to the wall. 
The difficulty with them is that the first grid point must be placed in the laminar 
sublayer (y'"~1) in order to properly connect the laminar sublayer with the log-law 
(wake) layer. Therefore, this kind of model requires the use of very fine grids 
near the wall boundaries. 

We selected the Low Reynolds number “kt-£t” model from Chien. For more 
information on this are see (Ref. 29) or (Ref. 30). The “kt-£t” equations have been 
modified to include the molecular viscosity that dominates in the near wall 
regions. The general form of the Low Reynolds number models is given by the 
following equations; 

|(pk) + ^(pu,k) = ^((p + + p(P - £ - D) (23) 

5t 5Xi dx, a. dX; 

J J J 


at 


(p£) + —(pUj£) = —((p + ^)—) + c J, ^ - C, f^ ^ + E (24) 
^ ’ dxr '' ax, ^ a/ax/ ' k ' k ^ ^ 


(25) 


The Chien model parameters appearing in the equations above are: 
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( 26 ) 


= 0.09, = 1.35, = 1.8, Qk = 1.0, a, = 1.3 

=1-exp(-0.0115y"),fi =1.0 
Rp 

f2=1-0.22exp(-(^)^) 

6 

OwU c 

D = ^,E = -2v(4)exp(-0.5y^) 

y y 

The turbulent kinetic energy was expressed using the turbulent intensity 
variation, which is provided by the following correlation (Ref. 29): 

l = 0.16Re ''® (27) 

for the core of a fully developed duct flow. 

As suggested in (Ref. 29), within internal flow problems, the duct hydraulic 
diameter is a good choice for a turbulent length scale. Further 
assumptions/settings included a) zero roughness conditions were set, b) to avoid 
numerical instabilities, the upwind first order scheme was used for advection, and 
c) the number of iterations of the solution process was increased in order to 
achieve the same convergence criteria as in the laminar cases. 

2. Transition Region 

As the Reynolds numbers increased, it was anticipated for the turbulence 
parameters to increase. Thus in each simulation, the variation of the at the 
wall boundary was examined. In order for y'' to remain in the order of unity further 
grid refinement was required. A y'' < 1 was achieved at the lower endwall 
boundary, but at the pins’ surface, y'' increased with height and with downstream 
location from the stagnation region as shown in Fig. 15. The larger values of y"" at 
increasing height were due to the extrusion power law that was used for the Z 
dimension during the grid construction. The increase in y"" along the pins’ 
perimeter can be explained by the acceleration of the flow due to the diffuser’s 
shape region between pins, thus causing a local increase of the turbulence 
levels. Figure 15 results were obtained at a Reynolds number of 3,500. 
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Figure 15. y"" variation with the half pin perimeter at several heights 


Despite further grid refinements and increased number of iterations, it was 
not possible to reduce the y"" values. Hence, turbulent models require very fine 
grids, a lot of computer memory, and even greater amount of computing time. 
Figure 16 represents a three dimensional carpet plot of y"" for a z-axis cut of the 
model for a Reynolds number of 3,500 at z=0. 



Figure 16. y"" variation in the lowest cells layer of the model 
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Also very encouraging from a heat transfer point of view is that the 
maximum heat fluxes at the pins’ surface are met at the 35-40 degrees 
downstream of the stagnation point. The in this region is still increasing, 
because it is the initial region of the flow acceleration. As expected from theory, 
very similar to the y"" trends were met for the turbulent kinetic energy. 

We try to identify the transition in a plot that compares the effective array 
heat transfer coefficient with the required specific fluid friction power. This plot fits 
our purpose well, since it combines both heat transfer and fluid flow 
characteristics. Figure 17 provides this information for both the laminar and the 
turbulent models for the Reynolds numbers between 1,000 and 8,200. 
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Figure 17. Comparison of laminar and turbulent models for FIX-9 


Both models overlap until Re > 2,000. After that point the laminar model 
under-predicts hgff. In the configuration of Flamilton (Ref. 7), the fall of the laminar 

model occurred at Re > 1,000. Therefore, first estimates indicate that transition 
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occurs in different regions for different configurations and depends on the area 
density of every array. For heat exchangers with smaller area magnification 
factors the transition region is shifted toward higher Reynolds numbers. At even 
higher Reynolds numbers shock waves are formed. This will be examined more 
extensively in the following chapter. 

The Nusselt number follows similar trends as the heat transfer 
coefficients, therefore Nusselt numbers decay near and after the shock regions. 
The shock wave problem was somewhat expected due to the micro-scale 
character of the exchangers. When comparing micro- with macro-scale heat 
exchangers at similar Re numbers, the micro-scale heat exchangers have to 
operate at higher velocities, which give rise to the shock problems. But since 
CFD-ACE is not a robust compressible solver, it is not very reliable for analyzing 
problems at these high velocities. Consequently, we use the results of the 
turbulent model until a Re < 6,200, where the maximum Mach number was 
approximately 0.83. 

Another aspect that has to be examined is the variation of the effective 
viscosity. We found that the effective viscosity reached and exceeded values five 
times greater than the laminar viscosity in the core of the flow (upper symmetry 
plane boundary condition) in a regime of the Reynolds numbers around 2,000. 
Figure 18 presents a carpet plot with the effective viscosity values for a z-axis cut 
of the array just below the upper symmetry boundary of the heat exchanger for 
the Reynolds number around 2,000. 
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Figure 18. Variation of effective viscosity in the HX-9 for Re=2000 

It is generally suggested that when this criterion has been reached, the 
modeling must switch into the turbulent regime. 

D. EXPERIMENTAL CORROBORATION 

1. Validation with Historical Data 

As discussed above, the current CFD results were compared to existing 
experimental data. Such experiments used the same configuration 
(X/D=S/D=2.5, H/D=1) from, Chyu (Ref. 10) and Metzger et al. (Ref. 6). Chyu et 
al (Ref. 12) conducted an experimental work based on a well known mass 
transfer technique, the naphthalene sublimation experiments, and through the 
use of the mass transfer analogy, obtained heat transfer results. The difference 
between Chyu and the current study are the following definitions: Chyu’s 
Reynolds number was defined as follows: 



where Vmax is the average flow velocity in the minimum flow area (between 
a pins row) and D is the actual pin diameter. 

Chyu defined the Nusselt number as follows: 
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where, again, D is the pin diameter and h is the heat transfer coefficient 
that comes from; 


h = q"(T.a„-T,,J(30) 

where q"is the heat flux and Tbuik is the bulk temperature in the channel 
based on the mass transfer analogy. Because it was impossible to equate their 
Tbuik using the information provided in their paper, and because the bulk outlet 
temperature used in the current study (the net outlet temperature of the coolant 
fluid,) was not consistent with theirs, it is better to compare their Nu number 
results with a slightly modified version of our Nu number (where the pin diameter 
is used instead of the hydraulic diameter.) The results can be compared by 
examining Table 7 and Figure 19. 


Table 7. Data comparison of present study with Chyu et al (1999). 


Reoh 

Nuoh 

1033 

19.7 

1498 

24.7 

2067 

30.0 

2790 

36.4 

3617 

42.7 

4650 

48.7 

6200 

51.8 

6613 

51.3 

8267 

36.9 


Present study Re 
equivalent 
969 
1406 
1939 
2617 
3393 
4362 
5816 
6204 
7754 


Present study Nu 
equivalent 

14.3 

17.9 
21.8 

26.4 

30.9 
35.3 

37.5 
37.2 
26.7 


Nu from Chyu et 
al correlation 

17.6 

21.9 

26.4 

31.5 

36.6 
42.4 
50.1 
52.0 
59.3 


The correlation that Chyu provides is Nu = 0.32Re°®®^ for Re>6000. The 
above extracted data is from an extension of their correlation toward lower Re 
numbers because our micro heat exchanger could not reach their higher 
Reynolds due to the shock wave issues. 
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Nu vs Re comparison with data from Chyu et a) (1999) 


■ Present study Nu equivalent | 
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Figure 19. Nu vs. Re comparison with historical data for HX-9 


If we consider that Chyu tested an array with 7 rows of pins and the 
current CFD model for FIX-9 had only 6, the difference in the definitions of the Nu 
number and other parameters (such as the viscous dissipation not being included 
in this study,) the agreement and the scale of the findings is sufficiently 
satisfactory. Moreover, the parallel trend that the two curves of the compared 
data are following is a very positive result. 

By comparing our findings with Metzger et al. we again achieved a very 
similar trend, with smaller differences of about 15%. Metzger used the same 
definition as Chyu et al. for the Re and Nu numbers, while his heat transfer 
coefficient was derived from the equation 


h = 


q/A 


wetted 


(^wall ■ ^ref ) 


(31) 


where for the W definition and computation the reader is referred to (Ref. 
6). The fact that the computational model’s trend-line develops in the region 
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between the two predictions of the two empirical researchers, is very positive for 
the performance fidelity of the model. 


Nu vs Re comparison witti Metzger [1982] data 
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Figure 20. Comparison of Nu as a function of Re with Metzger’s results 


Also Metzger in his study, found that the average Nusselt number 
decrease after the sixth row is negligible for this configuration and it can be taken 
as if it has reached a steady value of an infinite length array. 

In an earlier study Chyu (Ref 28), also researched the pressure losses 
and the friction factor for the same specific configuration. He defined the friction 
factor as follows: 


2AP 
pV ^ 

1^ max 


(32) 


where AP is the pressure difference between the inlet and exit sections of 
the array, N is the number of pin rows and the Vmax is defined as above. His 
findings are presented in figure 20. 
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Figure 21. Friction factor vs Re from Chyu [1989] for X/D=S/D=2.5, Fl/D=1 


The present study produced equivalent parameters depicted in figure 21 - 
agreement is excellent! Chyu’s experimental array consisted of 7 rows of pins 
and appears that this difference did not have the influence that the Nu number 
results showed. 


Present study f variation based on Chyus definitions 



Friction factor variation of the turbulent model for FIX-9 
43 


Figure 22. 











Marques et al. (Ref. 16) tested the same configuration as Chyu et al., but 
in the micro-scale. Marques’ pin diameter was 500pm, the same as to the current 
study. Accordingly, the HX-9 design is identical with the design of Marques, but 
with one row less. Also because Marques et al. heated only the one endwall, we 
could not compare the present study’s heat transfer findings with theirs. However 
a comparison of the friction factor could be very useful. After they corrected the 
parameters, they validated their results against Chyu’s study, as shown in figure 
23. 



Figure 23. Friction factor as a function of Re number from Marques [2004] 


Consequently, there is a very similar trend for the values of the friction 
factor predicted above which seems to tend towards a constant value at Re > 
6,000. This corresponds to a value of Re > 5,000 in the current study. 

2. Validation with NFS Experimentai Research 
Roussakies (Ref. 17) and Bosserman (Ref. 18) conducted experimental 
research using identical configurations to the current study. Roussakies worked 
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with configurations exactly similar to HX-1, 2, 10 while the Bosserman tested all 
of the current configurations except for HX-7,8 and 9 due to problems with 
blockage effects of the arrays’ passage ways. Therefore, only HX-1 and 10 were 
validated against their results. 

Bosserman’s results of the friction factor for some of the HXs like 3, 4 and 
6, in which he found numbers of 15 < f < 100 for Re=100, indicate trends that can 
only be explained by the heat exchangers being clogged (not properly released.) 

Figure 24 plots the Nu for the first heat exchanger for a range of Re. 



Comparison of Nu vs Re for HX-1 
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Figure 24. Comparison of Nu vs. Re for HX-1 with experimental work 


From the plotted data, we have good agreement for Re > 600, where the 
points are shown converging in a similar trend. We noticed a similar deviation 
with the Nu results for the average array heat transfer coefficient when it was 
plotted as a function of specific fluid friction power. This was expected, since both 
parameters are in great dependence. 
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Figure 25. Comparison of harray.ave vs E for HX-1 with experimental work 


Regarding the friction factor, our findings were again corroborated. Our 
results were found to be between the trends of the experimental results. Even 
more positive is the fact that the data showed similar decay as the Re increased 
- shown in figure 26. 


Comparison of f vs Re for HX-1 


26 

25 


1 B 


s IS 

t 

1 3 
^ 0 
06 
05 
03 
00 


I » Present compuatioftal model | 
I ■ Roussakies's results [2004) | 

I A Bou6rrt>an’Bf6iullft[2005| ] 


—I-1 •-1---1-'-1 • I-* 

100 200 300 400 500 000 700 600 000 1000 1100 


Re 


Figure 26. Comparison of f vs Re for FIX-1 with experimental work 
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The comparison of HX-10 with experimental data led to stronger 
corroboration, especially with Roussakies’s findings. His Nu variation with Re is 
closer to the computational model results than Bosserman’s. 


Comparison of Nu vs Re for HX-10 
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Figure 27. Comparison of Nu vs Re for HX-10 with empirical data 


Roussakies predicted values for Re > 400 close to the numerical model’s 
trend-line. In his data, the same behavior is seen in the average heat transfer 
coefficient values. This happens according to figure 28, for low friction power 
values, as compared to Bosserman’s results, which also predict a well 
established agreement for the friction factor. The trends between Roussakies 
and our model in figure 29 are extremely similar and the values are very close. 
The predicted friction factor values are almost identical for a region of Re > 400. 
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Figure 28. 
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Figure 29. Comparison of f vs Re for FIX-10 with empirical data 
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Consequently, the computational models of the heat exchangers that have 
been developed and are discussed in this chapter, produced results very similar 
to earlier experimental data in the macro-scale, where a suitable correction was 
applied, and extremely similar where the accurate equivalency from their 
definitions was known and has been implemented. As to the comparison with the 
micro-scale results, we are sufficiently satisfied with the existent corroboration. 

Therefore, based on the validation by both the micro and the macro scale 
data, the results of the macro-scale experiments can be directly applied to micro¬ 
scale heat exchangers. 
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IV. RESULTS-DISCUSSION 


Here we discuss selected performance characteristics of the heat 
exchangers that were analyzed (shown in Table 3) which for convenience are 
divided into two sets according to their pin diameter. The characteristics are 
chosen to highlight various features of the behavior of the heat exchangers. 

A. HEAT TRANSFER FLUX, RATE, AND COEFFICIENT 

The heat transfer coefficient is one of the most critical and interesting of 
the examined parameters. Since it was not possible to directly examine the 
behavior of the local heat transfer coefficient in CFD - VIEW, we inspected heat 
fluxes and heat transfer rates instead. 


Pins and enwalls performance comparison for HX-3 



Figure 30. Pins and endwall heat transfer comparison for HX-3 


HX-3 was found to be one of the best performers, since it had a very 
packed configuration but a bigger pin diameter. The contribution of the endwall 
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and the pins to the heat transfer, along the total length of the array is compared 
in figure 30 for Re=500. It can be seen that for the first row of pins the heat 
transfer rate achieved values more than 250 percent greater than the endwall, 
which was reduced greatly for the following rows without any great alteration in 
the percentage difference contribution. To confirm that result, we should inspect 
the heat flux distribution at those wall boundaries. The endwall showed very high 
fluxes for the very initial part of the Hx, but reduced very quickly and with a very 
high gradient along the stream-wise distance as is presented at the figure 31. 



- W/n^2 


Figure 31. Endwall’s heat flux distribution in HX-3 


After the first row the values of the actual wall flux were falling to a level of 
10,000 W/m^. On the other hand when we examined the heat flux according to 
pin chord-length in the first row for several heights, we noticed that it reached 
near maximum values from almost a quarter of the half height which was 
modeled. This means that the pin keeps the highest achieved heat flux values for 
more than 75 percent of its length. Furthermore, the maximum values for the 
heat flux reached are around 35-40 degrees downstream from the flow direction 
angle and are maintained until 100 to 105 degrees, from which point they start to 
reduce dramatically due to generation of vortex shedding and recirculation 
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effects. Figure 32 is a plot of this behavior with angle from the forward flow 
direction to 360 degrees around a pin. 


Heat flux distribution across pin perimeter for Isr row’s pin in HX-3 



Figure 32. First row’s pin heat flux distribution for several heights in FIX-3 


It is worth noting the scale of the heat fluxes presented and the great 
potential of the micro heat exchangers to dissipate large heat fluxes. The values 
of heat flux that the pins reach are extremely high and were highest for the 
second row as illustrated in figure 33. 
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Heat flux z-axis pfofile for hlX-3 in 35 downstream 



Figure 33. Heat flux z-axis profile for 35° downstream the pin surface 


After the second row, the values of the heat flux were reduced, but not 
with rates as high as those for the endwall. The second row was the one that 
clearly achieved the greater values and this is probably happen due to the 
acceleration of the flow after passing the first row, due to the nozzle effect which 
drives it more toward the next row’s pin center. In figure 34, which presents the 
contour distribution of the heat flux in a first and a second row pins, the red 
region of higher flux values has shifted toward the stagnation point for the second 
row pin compared to the one of the first row. However, the first row finally 
produced bigger heat transfer, when we look back at the comparison plot for the 
pins and endwall. This could be explained by the reason that the first row pin 
achieved greater area average heat flux values. 
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Figure 34. Heat flux contour plot of 1®* and 2'^'^ row pins in HX-3 


By averaging the heat flux values at these two pin surfaces using the 
calculator tool of CFD - VIEW we confirmed that result, but as we checked and 
for other arrays, we found that it did not happen in all the configurations. Also, 
because the great difference in the heat transfer participation, that pins and 
endwall resulted it was important to inspect how other configurations behaved. In 
a comparison for a Reynolds number equal to 500 for HX-1, we got opposite 
results to those in HX-3. In HX-1, the endwall provided more heat transfer, but 
not with as much difference as the pins did in HX-3. 
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Pins and endwalls performance comparison for HX_1 



Figure 35. Pins and endwall heat transfer comparison for HX-1 


The greater differentiation is found in the contribution of the first row, and 
as we checked the following rows, we noticed that this difference dies quickly, 
mostly due to the more general wall behavior. The pins continue to contribute a 
high enough amount, as compared to their initial amount of heat transfer. At this 
point, we should mention that HX-1 did not return the greater heat transfer 
coefficient, as compared to the other arrays. To get a total picture of this matter, 
we examined HX-2 which has the increased span-wise ratio S/D. The 
participation in this array was comparable for both - the endwall produced greater 
heat transfer in the first row only while the pins prevailed in the following rows 
with their more constant performance. 
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Pins and endwalls performance companson for HX-2 



Figure 36. Pins and endwall performance comparison in HX-2 


The pins produced a 20 - 30 percent greater heat transfer from row 2 to 
row 10. Consequently, the pins and endwall contribution is not constant; in fact, it 
is strongly dependent upon the configuration. More specifically if we examine the 
area magnification factor of the inspected arrays, we can state with great 
assurance that the pins’ participation depends on the area density of every array. 
Thus, in HX-3, with the higher area density, the pins contribute the greatest in a 
comparison of the three arrays. 

Furthermore, in both FIX-1 and FIX-2, the second row of pins produced the 
greater heat transfer rate in contrast to FIX-3. Accordingly, which row produces 
the maximum average flux and subsequently which row dissipates more heat is 
something varied, and depends upon the configuration too. Flowever, the heat 
flux in all of them reached its maximum local values in the second row. The result 
that the second row achieved the higher local heat fluxes is supported by figure 
37 which is a plot of the flux distribution at the upper symmetry plane along the 
pin chord-length for all rows of FIX-1. 
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Comparison of row's pin heat flux distribution for HX-1 
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Figure 37. Heat flux distribution for the pins of HX-1 at the core of the flow 


However, in order to make a more fair judgment of the pins and endwall 
contribution, we use a chart that illustrates the attained average heat fluxes of 
every row, which are independent of the flow wetted area factor, as in figure 38. 
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Figure 38. Heat flux comparison of pins and endwall in HXs1, 2, 3 
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Here we can observe that the difference in the first three configurations 
varies from 50 to 250 percent, based on the endwall values. Moreover, it shows 
more clearly that the performance of the pins is more generally dependent on the 
heat exchanger’s wetted area magnification factor. In addition, it is clear that the 
second row does not always accomplish the highest flux as happens in HX-3. 

For the average heat transfer coefficient for the two sets of arrays, figures 
39 and 40 present their variation with Reynolds number. It is clear that HX-3 and 
7 achieve the greatest values in the higher Re number regime for the two sets. 



Figure 39. Average heat transfer coefficient vs Re for HXs 1-4,9 


From a more total view, exchangers 6, 7, and 8 reached the higher 
coefficients, which might have been expected since they have the greatest area 
density. When we inspected the first set, HX-4 (X/D=3, S/D=1.25) was found to 
have greater performance, compared to HX-2 (X/D=1.25, S/D=3), suggesting 
that the stream-wise spacing variation proved to have a greater effect compared 
to span-wise distance. This result is consistent with the theory, since in HX-4, the 
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flow is forced to follow a more circuitous path around the pins and, as a result, to 
increase the interaction with them. In HX-2 with its wider passages, the gain in 
the heat transfer is not so significant, since the flow has greater direct through- 
flow regions at the side of the pins. Also worth noting are the values for HX-3 
which are far higher in the figure than all the others. The combined reduction in 
both spacings, stream-wise and span-wise, was the main thing that forced the 
third configuration to perform so much better than the others and to attain a 
similar performance as HX-6 and 8 with the small pin diameter. 
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Figure 40. Average heat transfer coefficient vs Re for HXs 5-8,10 


In the second set of Hx, the pin diameter was around three times less, 
which resulted in an increased H/D ratio. This was due to the constant height 
channel that we kept for all the configurations. In comparing the two plots of the 
two sets, it is easily noticeable that the main trends of the performance that have 
been commented on previously are qualitatively the same for the second set, but 
with a great change in between the configurations’ quantitative relation. This is 
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very significant, since it proves that the alteration of one of the three ratios, the 
pin height to diameter ratio, had a significant effect on the resultant performance. 
Thus the optimization of the micro heat exchangers is not so simple since it 
involves a complicated interaction of the different length ratios in the Hx. 

It is also important to note that the increase in the H/D ratio had an 
amplification effect on all of them, but with a different factor for every one of the 
exchangers as seen from the above discussion. 

A chart plotting the heat transfer coefficient versus specific fluid friction 
power, based on wetted area would be very useful in order to evaluate the 
performance. 


vs E for HXs 1-4,9 



Figure 41. Performance comparison for HXs 1-4,9 


The optimum Hx would be one that maximizes the heat transfer rate while 
minimizing the frictional losses in the flow. For the first set of arrays figure 41 
provides that information. Instead of comparing the average heat transfer 
coefficient, we selected the effective one, because it better represents the actual 
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heat transfer gain from every configuration, since it was derived by taking into 
account the magnification due to the array area density. The plotted points on 
each line are for increasing Reynolds number values of 100, 200, 300, 500, 700, 
and 1,000 in the direction of increasing E. 

Between HXs 1, 2 and 9 there is no great difference, except that the 
second one which, for the same Re values, reaches slightly higher values of heat 
transfer coefficient, with the analogous required energy, than the others would 
demand operating at a similar point. However, HX-3 and HX-4 have a great 
difference between them, with the fourth starting from a value for the hgff of 450 
W/m^/K and reaching a maximum of 1,000 and the third starting from that point 
and reaching the amount of 4,000 W/m^/K for the same Reynolds number. If we 
consider that, for the heat transfer coefficient near 1,000, one demands 100 
times more power to drive the flow, we can see that there is a great difference 
between their performance. 

HXs 2 and 4 are closer with the fourth one performing better until a value 
of Re equal to 500, whereas after this point, HX-2 produces slightly (maximum 30 
percent) smaller values of heat transfer coefficient while demanding around three 
to four times less energy to operate. 

In the second set of inspected configurations, all the plotted performances 
have moved toward the right upper region of the chart except for HX-10. 
Furthermore, it is interesting that their differentiations are scaled similarly in both 
figure’s directions for all of them except HX-7. The most impressive characteristic 
of this comparison is the heat transfer coefficient values that HX-7 reached 
considering its higher area density. But the tremendous amount of fluid friction 
power needs to be compared to all the other arrays. The deterioration in its 
performance was clearly noticeable when it approached the shock wave 
operation region, as the heat transfer reduced while the required specific power 
grew from 200,000 to almost 500,000 W/m^. 
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Figure 42. Performance comparison for HXs 5-8,10 


To increase the heat transfer performance, it is much more cost effective 
to move upward in the diagram to a more packed configuration, instead of 
changing the flow attributes and increasing the Reynolds number. 

Also, if we consider HX-1 and 10, which have almost identical 
performances but differ a lot in their geometrical characteristics, the gain from 
increasing the H/D ratio in the tenth has been diminished by the increase of the 
other two spacings, X/D and S/D, as compared to the first array. This is an 
indication that the height to diameter relation has a stronger influence compared 
to the other two ratios, and thus needs more attention in the optimization 
process. 

Another way to inspect the heat transfer performance of the tested 
configurations is to compare them against the performance of the similar 
dimensions empty duct. Thus in this case, it can be examined the difference that 
the presence of the pin made, not only by contributing with their surface but 
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enforcing the flow to follow the more tortuous paths around them and increase 
the heat dissipation of the wall boundaries. 
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Figure 43. Variation of arrays’ effective heat transfer coefficient over the empty 

duct for HXs 1-4,9 


Figures 43 and 44 provide the ratio of the effective array heat transfer 
coefficient over the corresponding of the empty duct for various Reynolds 
numbers. The empty duct heat transfer coefficient came up after taking in 
account the constant Nusselt number of the infinite width square duct, which is 
equal to 7.54 for the thermally developed flow, as it is reported by Kays and 
Crawford. 

The configurations like HX-3, 6, 8 achieved from five to almost thirty times 
greater heat transfer coefficients when consider their corresponding empty duct, 
but the HX-7 returned values of the order of thirty to more or less one hundred. 
This, as we will examine the friction factor in a following section, was not without 
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an admirable energy requirement that almost forbids any thoughts of applicability 
that can this configuration would provide. 

Furthermore, the fact that for the arrays the heat transfer coefficient 
increases significantly by increasing the flow characteristics (Reynolds number), 
which does not occur so intensively for the case of the empty duct, is also 
important. 



Figure 44. Variation of arrays’ effective heat transfer coefficient over the empty 

duct for FIXs 5-8,10 


B. NUSSELT NUMBER VARIATION 

The variation of Nusselt number with the Reynolds number is shown in 
figure 45 for the first set of arrays. The trends are not very different compared to 
the behavior of the average heat transfer coefficient. 
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Nu vs Re loi HXs 1-4,9 



Figure 45. Nusselt number comparison for HXs 1-4,9 


The actual relation between the heat exchangers has changed a little 
because the Nusselt number has been influenced by the hydraulic diameter, as 
defined in a previous chapter, according to equation 2. This differentiation first 
concern HX-2, which has moved in a lower scale, and HX-9, which has moved 
upward relative to HX-1, always compared with the heat transfer coefficient 
trends. This is valid with the variation of their hydraulic diameter, which for HX-2 
with the smaller volume and the same cross sectional area as HX-1, results in a 
smaller diameter. For HX-9 the reduction of both ratios (X/D and S/D) produced a 
smaller fractional reduction in the volume than in the cross sectional area, and 
consequently, the enlargement of the hydraulic diameter. Furthermore, when 
considering the Nusselt number variation, the opening between the HX- 3 and 
HX-4 has been decreased, since the hydraulic diameter of 567pm of HX-4, 
compared to the 307pm of HX-3, appears as an amplification factor for HX-3. 
This trend is clearer in the second set, where the increase of the H/D ratio 
increased the heat transfer coefficient for both HXs 7 and 8, but the further 

increase to the hydraulic diameter for HX-8 (Dh,8=354pm vs Dh,7=137pm) 
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dominated and caused not only the gap’s total diminishing but also the eighth 
array’s outperformance over the seventh until an Re equal to 300. 

By continuing with the exchangers in the second set, as was done for the 
heat transfer coefficient, the performance curves in figure 46 include the effect of 
the increase in the H/D relation and the alteration to their trends by our 
characteristic length, the hydraulic diameter. 
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Figure 46. Nusselt number comparison for HXs 5-8,10 


In judging all the configurations’ performance considering the Nusselt 
number, against the heat transfer performance of an open channel with a 
constant surface temperature boundary which produces an average Nu equal to 
7.54 (Ref 26), almost all of them attained much higher values for the whole 
tested range of the Re number. The only exceptions were HX-1 and HX-2, which 
cross the limit of the constant open channel performance around a Reynolds 
number of 200. 
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C. FRICTION FACTOR EFFECTS 

The most obvious effect that occurred in the friction factor, in comparing 
the two sets, was the magnification due to the increase of the H/D fraction that 
resulted from the reduction of the pin diameter, while keeping the channel height 
constant. The reduction of the pin diameter caused a decrease in the hydraulic 
diameter and also a great increase in the average array velocity. However, these 
two changes were not enough to offset the decrease of total array length and, 
more significantly the increase in the pressure gradient between the entrance 
and the exit of the exchanger. 


Friction (actor (vs Re for HXs 1-4.9 



Figure 47. Effect of Reynolds number in friction factor of HXs 1-4,9 


Furthermore, for the friction factor, like the other parameters, we note that 
the amplification due to the alteration of the H/D was not proportional for all the 
arrays. Thus, this is another proof that the changes to the geometric fractions do 
not act independently and must be considered together. 
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If we consider the difference in the effects between HXs-2 and 4, or 6 and 
8, since they have the same hydraulic diameter, and thus the same average 
velocity, the main cause of it, is the great difference between the two differential 
pressures of the compared arrays. Moreover, this cause is consistent with the 
theory given that the smaller S/D ratio for HXs-4 and 8 produces a narrower 
span-wise spacing for the pins and consequently the flow streamlines negotiate a 
more serpentine path around the pins, which leads to the greater pressure drops. 
On the other hand in HXs-2 and 6 there’s a more direct path and more smooth 
streamlines, and result in smaller energy amounts being needed to drive the flow 
through their interior. 


Friclion factor f vs Re for HXs 5-8.10 



Figure 48. Effect of Reynolds number in friction factor of HXs 5-8,10 


Another aspect worth noting is the high friction factor that HXs 4 and 8 
produce, especially as compared with HXs-3 and 7. The main reason for that is 
the dependence on the squared of the average velocity at every Re number and, 
the increased hydraulic diameter. 
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By comparing the pressure drop across the inlet and the exit of every 
configuration, with the pressure drop required from the corresponding 
dimensions empty of pins duct, to run the flow through it, we get figures 49 and 
50. The amount of energy required for the several arrays in order to increase the 
heat transfer performance shows similar behavior to the corresponding friction 
factor. 



Figure 49. Ratio of HXs 1-4,9 pressure drop over the empty duct case vs. Re 


When considering an optimization and a selection process of the most 
suitable HX for an application, the amount of the pressure drop required for 
returning the maximum heat transfer performance or the inverse, has to be 
examined. Thus, HX-7 requires around five times more pressure difference of the 
increase in heat transfer that provides, HX-8 proved to require six times more, 
while HX-6 proved to be more economical by demanding three times greater 
pressure difference than the times that increase the heat transfer when 
compared to the open duct case. 
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Figure 50. Ratio of HXs 5-8,10 pressure drop over the empty duct case vs. Re 


D. MACH NUMBER VARIATION 

The deterioration of the heat transfer performance of HX-7 is a very 
interesting manifestation of Mach number variation. If we also consider the 
extremely small values of the hydraulic diameter for these micro-scale devices it 
is obvious that high velocities have to be simulated to provide the desired tested 
Reynolds numbers. 

Recognizing that CFD - ACE is not the most suitable solver for studies 
that involve Mach numbers over 0.7-0.8 the presented results and effects are 
rough estimates and must be seen as only adequately accurate in revealing the 
main trends. A compressible solver like CFD - FASTRAN would be required to 
resolve shock behavior with greater accuracy but was not used here since it was 
peripheral to the main study. 

However by using the CFD - VIEW, we tried to identify the maximum 
Mach numbers that occurred in every array and to plot them with respect to the 
Reynolds number. They are presented in figures 51 and 52. Suitable correlations 


71 












to define this operational region would be useful in the design process. The 
coefficients of a power law correlation below are shown in table 8. 


Figure 51. 


M = C,Re^^ (33) 


Maximum Mach vs Re for HXs 1-4.9 



Maximum Mach number variation with Re number for HXs 1-4,9 


Maximum Mach vs Re for HXs 5-8.10 



Maximum Mach number variation with Re number for HXs 5-8,10 
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Figure 52. 


























We have to mention again that these correlations are not intended as an 
accurate prediction of the flow Mach number, but rather as a rough guess in 
order to assist in an attempt to avoid any interference with shock waves. 


Table 8. Mach number correlations’ coefficients 


HX # 

Ci 

C2 

1 

0.0002 

0.9644 

2 

0.0002 

1.009 

3 

0.0006 

0.9445 

4 

0.0003 

1.0442 

5 

0.0002 

0.9635 

6 

0.0003 

1.0158 

7 

0.0005 

1.0179 

8 

0.0011 

0.991 

9 

0.0002 

0.968 

10 

0.0002 

0.9606 


HX-7 demonstrated these shock conditions at the high inlet velocities that 
the theory dictated for the tested Re numbers and found to have maximum Mach 
number values of about 1.1. Furthermore, if we consider the inlet velocity of 
60m/sec, or the average of 121 m/sec, or the maximum value of the fully 
developed profile of 90m/sec for the corresponding Re equal to 1,000, and the 
drastic reduction of the available flow area between the pin passageways, shock 
waves in these nozzles like sections seem unavoidable. The contours of the 
Mach number and pressure distribution around the throat, illustrated in figure 53, 
is another indication. We recognize that the coarser unstructured grid around the 
centerline is not appropriate for capturing the occurring shock conditions. 
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Figure 53. Contour plot of pressure and Mach number distribution in the 4*^ 

row of HX-7 


In addition, the fact that the critical pressure ratio of 0.528 when a shock 
can occur is reached is a clear indication. For all these reasons, we are confident 
that the shock behavior was captured to some degree by ACE and is reflected in 
a decline in the heat transfer coefficient and in Nusselt number values around a 
Reynolds number of 1,000. 

Bosserman (Ref. 18) observed a very similar condition for FIX-3, the 
denser of the first set in his experimental research. Also FIX-3 contains the same 
area percentage reduction between the pins as FIX-7 does. Fie observed similar 
decreased performance conditions around a Reynolds number of 1,000 at the 
same parameters. Figure 54 presents an indication of the variation of Nusselt 
versus Reynolds number for FIX-3 with his data. 
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Figure 54. Nu with Re variation for HX-3 [from Bosserman, 2005] 


In considering the data in table 9 that provides the required velocities, inlet 
duct, and average velocities, we concluded that a shock condition must have 
happened at a smaller Reynolds number than expected, since at a Re equal to 
700, the theory gives an inlet velocity of 17m/sec for which it is very difficult for 
the flow to be chocked under the geometrical features of this array. However, if 
we consider the blockage affects that Bosserman faced as well, it is very 
possible that the flow was at higher velocities and that he experienced shock 
conditions that are consistent with the gathered data. 


Table 9. Inlet duct and average array velocities for HX-3 for several Re 


RSoh.ar 

m dot ar(kg/s) 

Uave ar(m/s) 

Vin(m/s) 

100 

1.5463E-06 

5.39 

2.68 

200 

3.0926E-06 

10.78 

5.36 

300 

4.6389E-06 

16.18 

8.05 

500 

7.7315E-06 

26.96 

13.41 

700 

1.0824E-05 

37.75 

18.77 

1000 

1.5463E-05 

53.92 

26.82 
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E. SLUG VELOCITY PROFILE CONSIDERATIONS 

One of our considerations was whether the performance of the heat 
exchangers would change and if so how much, by applying a more realistic and 
not fully developed velocity profile as an input boundary. Therefore, the problem 
changed to a thermal constant surface entry length combined with a flow 
developing condition. Also it is well known (Ref. 26) that the thermal layer 
generally develops faster than the velocity boundary layer for fluids with a Pr < 5. 

Thus, we chose five arrays from both groups, with the smaller and the 
bigger pin diameter, in order to check the effect of all the spacing ratios, testing 
them in the new conditions. We removed the entry duct from the model and we 
applied directly a “slug” constant velocity at the inlet of the exchanger that 
corresponded to the desired Reynolds number. It is very clear from figure 55 that 
the performance did not change too much when the slug inlet velocity models 
were compared with the fully developed velocity profile models. We can observe 
a slight increase in the heat transfer coefficient at all the ranges of the Reynolds 
number, with another slight increase for the energy required to drive the flow 
through the exchanger. 
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vs E for comparison of HXs 2-4,6.8 


Figure 55. 
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The resulting behavior is consistent with the fluid physics, since from the 
flow perspective at the initial part of the endwall, the friction gradient will be very 
high until the flow will be developed relative to the endwall, resulting in increased 
fluid friction power and also in an enhanced heat transfer performance. 

As far as the Nusselt number is concerned, it is higher, especially for the 
lower x+ values. Also, as the x+ increases, the Nusselt number has to converge 
to the values of the thermal length only problem, since the effects of the 
combination of the two developing conditions are diminished. Thermal 
development is attained when the x+ coordinate reaches values of 0.1. 


Table 10. x+ values for HX 2-4,6,8 in Re=100,1000 


HX# Dh ar(m) L(m) 

2 0.000567 0.00625 

3 0.000307 0.00250 

4 0.000567 0.00600 


x+forRe=100 x+forRe=1000 
0.312 0.031 

0.230 0.023 

0.299 0.030 
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6 

8 


0.000354 0.00208 

0.000354 0.00200 


0.166 

0.160 


0.017 

0.016 


In table 10 we notice that for all the tested exchangers, the thermal 
development was reached for Re=100, however we are slightly below this point 
for the higher Reynolds number flow conditions. Thus, we expected the deviation 
in the Nusselt number between the similar geometrically models compared, to 
slightly increase as the Reynolds number was increasing, which did not occur. 


Nu vs Re for Comparison of HXs 2-4,6,8 



Figure 56. Nusselt number variation for compared velocity profiles in HXs 2- 

4,6,8 


The most probable explanation is that the performance characteristics had 
changed so little that it was of the order of the simulation accuracy, and we could 
not identify the above mentioned characteristic behavior. Also we noticed that the 
selection of the total number of rows for modeling every configuration was 
consistent with the corresponding x+ values in order to reach a near thermal 
developed condition, as was explained in the chapter two. 
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A plot of friction factor with Re in figure 57 helps to compare the role of 
inlet velocity profile on pressure drop. It is clear that there is very little influence of 
the velocity profile and noticeable only for HX-2 and HX-6. 


f vs Re tor comparison at HXs 2-4.6,8 



Figure 57. Friction factor variation for compared velocity profiles in FIXs 2- 

4,6,8 


The most logical explanation is that for the heat exchangers that provide a 
larger cross sectional area for the flow, the effects of the change in the inlet 
velocity profile were more important due to greater interactions with the endwall. 
In the other configurations, probably the presence of the pins and their dominant 
“denser” arrangement did not allow great changes in the compared parameters. 
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V. CONCLUSIONS AND RECOMMENDATIONS 


The objectives of this study were successfully met by computationally 
simulating the performance of ten micro heat exchangers with varied X/D, S/D, 
and H/D ratios. The code CFD - ACE was found to be very robust and user- 
friendly and was used with great confidence for the subsequent computations of 
the laminar flow regime and the performance maps of the interesting parameters. 
Most of the significant conclusions have been pointed out during the discussion 
in each section/chapter. 

In summary, we can say that the laminar modeling was shown to be 
accurate until a Reynolds number of 1,000. At higher values of Re there’s a 
likelihood of transition to turbulence, although the transition value is not constant 
and is greatly dependent on the geometric configuration of each array. 

The ninth and tenth configuration was well enough validated by historical 
macro-scale and current experimental micro-scale data. Thus, we can conclude 
that direct application of the macro-scale dimensionless parameters can be 
performed in the micro-scale dimensional heat exchangers in order to examine 
the performance of any desired configuration. 

Moreover, we noticed that the effects of an alteration of any of the 
geometric ratios influence the resulting effects of the variation of the other ratios. 
This was more pronounced for the H/D ratio. 

The performance maps of the interesting parameters like the Nusselt 
number, average and effective heat transfer coefficient, specific fluid friction 
power, friction factor, and Mach number were also provided for in all the 
configurations. They can be very useful tools during a design process. The 
tremendous ability and potential for heat transfer removal, especially for the more 
packed arrays, has been pointed out and also the extreme fluid power that some 
of them required. 
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An examination of the inner part of the compact heat exchanger was 
conducted, and, mostly, the behavior of the pins compared to the endwall was 
performed, which proved the great contribution that the pins provide in the total 
heat transfer removal process. 

The desired tested Reynolds number simulations in the micro-scale 
dimensions required much higher velocities which led to shock wave conditions 
very early for HX-7, Thus, the computation with varied properties was proved, the 
best selection. Also the presence of vacuum effects at the exit of some 
exchangers further enhanced more the choice of the variable properties 
modeling. 

The alteration of the inlet velocity profile that was performed by a constant 
“slug” velocity profile, in order to inspect the change in the parameters’ trends, 
did not return any considerable deviation of the initial maps. 

For all these reasons, the optimization and the design of micro-scale 
compact heat exchangers is a complicated process. It must be executed with 
considerable attention, since one has to consider more physical phenomena than 
in macro-scale design. 

Likewise to get a more total picture of the performance in future work we 
recommended that the viscous dissipation effects be taken into account in a 
computational model, in order to check how this negatively affects the output. 

Also the examination of configurations with other shape fins that will 
provide better flow characteristics and will be more effective is emerging, while 
the difference on the heat transfer performance that the liquids will provide has to 
be considered. 

Finally of a lesser concern is the fact that in some of the exchangers 
tested (only two to three) for the higher Reynolds number (700 and more for 
1,000), the residuals had a periodic pattern. This feature is an indication of 
unsteadiness, which we recommend should be checked by running an unsteady 
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simulation to inspect any effects on the final performance behavior of the heat 
exchangers. 
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